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The fission of highly charged sodium clusters with fissilities X > 1 is studied by ab initio molecular 
dynamics. Na24 4+ is found to undergo predominantly sequential Na3 + emission on a time scale of 
1 ps, while Na24^ + (5 < Q < 8) undergoes multifragmentation on a time scale > 0.1 ps, with Na + 
increasingly the dominant fragment as Q increases. All singly-charged fragments Na„ + up to size 
n — 6 are observed. The observed fragment spectrum is, within statistical error, independent of the 
temperature T of the parent cluster for T < 1500 K. These findings are consistent with and explain 
recent trends observed experimentally. 
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Starting with the work of Sattler et al. |]J on van der 
Waal's clusters, the study of the fragmentation (fission) 
of charged atomic or molecular clusters has proved a valu- 
able experimental tool for investigating the intrinsic sta- 
bility and binding forces of these objects. These stud- 
ies and others yield important insights into the behavior 
of matter at the small size limit and the size-dependent 
evolution toward bulk properties. Much recent work on 
fission has been devoted to metallic clusters, both ex- 
perimentally |g-[ll]] and theoretically g,|,[L||L6|] . Fission 
of metallic clusters is particularly interesting on account 
of the similarities and differences with the nuclear fission 
process 12 . Over a century ago, Lord Rayleigh [jlTj stud- 
ied the problem of the time development of the Coulomb 
instability of a surface-charged liquid drop, and conjec- 
tured that the excess charge would be removed by the 
emission of jets, rather than by fission into two parts of 
more or less equal size (as occurs for nuclei, where the 
charge is distributed uniformly through the volume). 

According to the Rayleigh criterion Jjl]] , a charged liq- 
uid drop is unstable against Coulomb forces when its 
fissility X > 1, where X = Ecou\/ {2E SUT f ) is propor- 
tional to the ratio of the Coulomb to surface energy of 
the drop. For a metallic cluster Mn® + , assumed spheri- 
cal with radius rgN 1 / 3 (r s is the Wigner-Seitz radius for 
the metal), this gives X — (lQirr^cr)" 1 Q 2 / 'N , where a is 
the surface tension, or X « 2.5Q 2 /N for Na7v^ + - Now, 
fission experiments on metallic clusters where the cluster 
is charged by laser ionization |2|-|7j have so far produced 
only clusters with X < 1, for which an energy barrier ex- 
ists against fission. Thermally activated fission may be 
observed, however, and the dominant charged fragment is 
found to be Na3 + , which has a closed electronic shell and 
is particularly stable; Na + has not so far been observed. 
In a different type of experiment [p[-pl|, sodium clusters 
are ionized by collision with a beam of highly charged 
ions, a technique that is expected to allow study of a 
much wider range of X and cluster temperatures. Coinci- 
dence measurements reveal multifragmentation processes 



in some cases, often with Na + as the dominant fragment. 

In this Letter, we offer the first systematic dynamical 
study of metallic cluster fission in the regime X > 1 using 
ab initio molecular dynamics (MD) jl8). We find that for 
X close to unity, sequential emission of mainly Na3 + is 
the dominant decay channel, while for multifrag- 
mentation with Na + as the dominant fragment occurs. 
Our simulations provide detailed spatial and temporal 
information on the fission process, and explain of some 
of the trends observed experimentally. 

On each time step of the cluster dynamics, we compute 
the density of valence electrons (and hence the forces on 
the ions) within the temperature-dependent Kohn-Sham 
(KS) formalism p9| , using the local density approxima- 
tion (at zero temperature) for the exchange-correlation 
functional E xc [p\. We use a real-space finite-difference 
method, recently developed by us incorporating a 
novel system of adaptive simulation cells that surround, 
adapt to, and follow each distinct fragment during a mul- 
tifragmentation, permitting the simulation to continue 
efficiently up to large fragment separations. To achieve 
better fragmentation statistics, at the expense of some 
loss of first-principles accuracy, we employ a soft, phe- 
nomenological pseudopotential |2l[| , which permits a rel- 
atively large real-space grid step A = 1.35 a.u. We do not 
expect our choice of functional E xc [p] or pseudopotential 
to affect significantly the main results for barrierless fis- 
sion. 

To study the fragmentation of a single species Na^v^" 1 " 
at an "initial temperature" Ti n , we run M dynamical sim- 
ulations arising from an ensemble of M initial conditions 
constructed as follows: (i) We optimize the geometry of 
the neutral cluster Na^v at K; (ii) we perform an MD 
run of > 20 ps for Na^ at 400 K (which is roughly the 
temperature of the clusters Na^ output by the cluster- 
source in the collision experiments 0-[ll]]); (iii) we take 
M ionic (nuclear) configurations {Rj } with velocities 
{Vj"-*} (n — 1-M) at regular intervals from this simula- 
tion; (iv) for each {R^™'}, we remove Q electrons and 
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FIG. 1. Fission barriers for two fission channels of Naio 2+ , 
as a function of the separation of the centers of mass of the 
two fragments, with zero corresponding to the initial cluster. 

re-equilibrate the remaining electrons to an electronic 
temperature T c \ = T; n ; and (v) we start the dynamics 
with configuration {R}"'} and velocities {AV^}, with A 
chosen to give an ionic (kinetic) temperature Tj on = 2] n . 
While the fragmentation spectra may depend to some ex- 
tent on the initialization procedure, the above procedure 
has been chosen to approximate the heavy-ion collision 
experiments |q-|ll]]. The collision time is fast, of order 
10 fs, so the ionic coordinates {Rr} are effectively frozen 
during the ionization process, as above. The scaling fac- 
tor A in the initial conditions is intended to approximate 
the extra energy "injected" into the ionic system by re- 
laxation of the valence electrons, which are excited during 
the collision. Steps (iv) and (v) above effectively assume 
this relaxation to be very rapid. We discuss the relative 
sizes of relaxation and fragmentation times below. 

Although we are mainly concerned with barrierless fis- 
sion X > 1, to assess the accuracy of our KS approach, 
and to make contact with previous experiment and the- 
ory, we first consider briefly the binary fission of small 
doubly-charged clusters. Our lowest-energy geometries 
for Najv (4 < N < 8) agree with previous density- 
functional theory (DFT) studies |22|Jl3i| , and our pseu- 
dopotential is adjusted to give a good fit to ionization 
and atomization energies for this size range, to within 
the scatter of previous DFT results ||^,[r|. We find 
doubly-charged clusters Najv 2+ to be unstable for N < 6, 
in agreement with Ref. |l3f| . In dynamical simulations 
with M = 10 initial conditions for each N and with 
T in = 600 K, Najv 2+ (5 < N < 10) undergoes binary fis- 
sion with Na3 + as the dominant fragment, in agreement 
with experiment (2]|| and previous theory Fission 
products other than Na3 + are found in only two cases: 
Na6 2+ — > Na5 + +Na + with about a 20% branching ratio, 
and Naio 2+ — v 2Nas + with about a 40% branching ratio. 

To understand the result for Nain 2+ further, we show 
m Fig. the energy barriers for the two observed fission 
channels, obtained by constrained energy minimization 
with the separation of the centers of mass of the two 



fragments specified. A double-humped fission barrier is 
found for each channel, and in the dynamics the cluster 
elongates to a "precursor state" where it may remain for 
several ps before fissioning, as also found previously in 
Rcf. Our barrier height for Na3 + emission is about 
0.5 eV, in reasonable agreement with the 0.7 eV found 
in Ref. [ fL3[ . The discrepancy may be due in part to our 
phenomenological pseudopotential. 

Let us now turn to a systematic study of the effect 
of fissility on the fragmentation spectrum for fissilitics 
greater than one. We shall consider the fragmentation of 
Na 24 Q+ for Q = 4-8 (X = 1.7-6.7) with initial tempera- 
tures T in = 400 K, 800 K, and 1500 K. We find Na 24 Q+ 
to be unstable (at K) in our KS model for Q > 4; 
Na24 3+ (X = 0.94) is just stable, with a barrier of about 
0.2 eV for removal of Na3 + . For each Q and Ti n , we run 
M = 10 simulations, each lasting up to 5 ps (Q = 4), 
3 ps (Q — 5-6), or 2 ps (Q — 7-8). A distinct final- 
state fragment Na„ 9+ is considered to have formed when 
all n ions in it are separated from the remaining ions by 
more than a cutoff distance r cut = 14.0 a.u. Its charge 
q is calculated as the total charge inside a box centered 
on the fragment with a border of at least 7.0 a.u. from 
any ion. Usually, q is integral to better than a few per- 
cent, and the identification of the fragment is unambigu- 
ous. But this is not guaranteed by the KS formalism: 
when two virtual orbitals centered on different clusters 
are nearly degenerate and overlap, the resulting KS or- 
bital may "split" between the two centers yielding frac- 
tional charges. This turns out to be particularly likely 
to happen with monomers or dimers emitted toward the 
end of a multifragmentation process. 

In such cases, it is usually possible to assign integral 
charges unambiguously by interpreting the electron wave- 
functions statistically. A typical example would be 

Na 24 8+ - 4Na+ + Na 18 2 - 85 + + Na°- 53 + + Na°- 70 +. (1) 

The first four fragments emitted here are Na + with 
very close to integer charge, but when the simulation is 
stopped (here after 2 ps), the remaining fragments are 
fractionally charged. We round the large fragment up to 
Nai8 3+ , and assume that the total remaining charge of 
1.08 shared by the two monomers is to be interpreted, in 
a statistical sense, as Na + Na + , with the probability for 
finding the charge +1 on a particular monomer given by 
the fractional charges. In this way, we often find neutral 
monomers or dimers emitted in the final stages of a mul- 
tifragmentation (but never among the initial fragments), 
which we interpret as evaporation from a hot residual 
fragment. Note that the slight excess of positive charge 
J2i Qi = 8.08 on the r.h.s. of Eq. (0) is due to spillout 
of electron density from the boxes used to calculate the 
total charge, which results in a slight underestimate of 
the negative electronic charge. When simple rounding or 
charge redistribution among equivalent species does not 
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FIG. 2. Snapshots of the Coulomb fission of Na24^ + for 
Q — 4-7 and an initial temperature Tin = 800 K. The time t 
for each snapshot (with t = the initial condition) is shown. 
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give a clear assignment of integral charges, we discard 
the simulation, which was the case for less than 5% of 
simulations. 

Some typical snapshots of Coulomb fission processes 
illustrating the main points are shown in Fig. ||. For 
Q = 4 (X = 1.7) (and also for Q = 3, X = 0.94), we 
find mainly sequential emission of Na3 + on a 1 ps time 
scale, with only rarely Na + emission. Such emission con- 
tinues until the large residual fragment (which ultimately 
develops a fission barrier) is too cool to emit further frag- 
ments, at least on the time scale of our simulation, here 
5 ps. (It is possible that with a longer simulation time 
we would occasionally observe emission of an additional 
singly-charged fragment.) As Q increases, we find emis- 
sion of increasing quantities of Na + , which is the domi- 
nant fragment for Q > 5. All singly charged fragments 
up to size six are observed in some quantities. 

Figure || shows in detail a disintegration of Na24 7+ 
(X = 5.1) as a function of time. In Fig. ||(a), taken 
at t = 0.10 ps, two Na + at the top and bottom left are 
about to leave the cluster. After 0.2 ps [Fig. ^(b)], a 
third Na + and a Na3 + start to leave from the rear. At 
t = 0.3 ps [Fig. ||(c)], the residual cluster is already highly 
deformed, and eventually emits a further Na + and Na3 + 
after about 0.7 ps, leaving a Nai4 + residue that remains 
stable up to 2 ps, when the simulation was terminated. 
The first few Na + ions emitted at about t ~ 0.1 ps were 
initially at the surface, where the valence electron den- 
sity is low (the excess positive charge tends to be located 
near the surface in a metallic cluster). We conclude that 
these ions were sufficiently weakly bound that they they 
simply accelerated outwards starting at t = 0. 

The precise fragmentation pattern found in a given 
run for large X is very sensitive to the initial condition, 
even for a given T; n . However, averaging over initial con- 
ditions, we find the mean number of a given fragment 



FIG. 3. Coulomb fission of Na24 7+ at times of (a) 0.1 ps, 
(b) 0.2 ps, (c) 0.3 ps, and (d) 0.45 ps. Shown are isocontours 
of the valence electron density, with Na + ions in the interior. 



emitted to be, within statistical error, essentially inde- 
pendent of T in for T in = 400 K, 800 K, and 1500 K. 
This is perhaps not surprising, since the Coulomb en- 
ergy of Na24 C ' ?+ is EqouI ~ 1-2Q 2 eV, while the ionic ki- 
netic energy is only E^n ~ 0.3 eV per 100 K and is thus 
small compared to -EJcoul for all cases considered here. 
Therefore, in Fig. ^ we have combined our results for the 
three T- m , and show the average number of singly-charged 
fragments Na„ + (n = 1-5) emitted per fragmentation, 
averaged now over 30 initial conditions for each parent 
charge Q. The most striking trend observed in Fig. [J is 
the smooth increase in number of Na + emitted with in- 
creasing Q, accompanied by a reduction in the number of 
Na 3 + emitted. For Na 24 4+ (X = 1.7), Na 3 + dominates; 
for Na24 8+ (X = 6.7), on average 20 times more Na + are 
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FIG. 4. Average number of Na n + fragments emitted per 
event versus the size n of the fragment, for the Coulomb fission 
of Na 24 Q+ (Q = 4-8). 
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emitted than Na 3 + . 

The average time elapsed before the emission of the 
first fragment varies from 0.1-0.2 ps for Q = 6-8, where 
the first fragment is nearly always Na + , to 0.4-0.6 ps for 
Q = 4-5, where Na 3 + is most often the first fragment. A 
combined electronic-ionic dynamical study within time- 
dependent KS [|l6| suggests that electron-ion relaxation 
times r r may be of order r r < 100 fs, and thus competi- 
tive with the time scale r/ of some faster fragmentation 
processes (17 > 50 fs). However, as we have seen, these 
faster processes involve immediate acceleration of Na + 
away from the surface region of the cluster, a process 
which does not require electron-ion relaxation to occur. 
Moreover, our mean fragment spectra are found to be es- 
sentially independent of Ti n up to at least T; n = 1500 K. 
We do not believe, therefore, that our approximate treat- 
ment of electron-ion relaxation via the parameters A and 
Tj n will lead to a significant qualitative error in our re- 
sults. 

The trend shown in Fig. [| is consistent with the obser- 
vation in a heavy-ion collision experiment jl(]] that Na + 
dominates the inclusive small fragment spectrum when 
Xe 20+ is the projectile, while Na 3 + dominates when Ar 3+ 
is the projectile, since the heavier ion Xe 20+ should pro- 
duce clusters NaAr^ + with higher charge Q. Our results 
are consistent also with another collision experiment [ [il"| , 
in which Na 3 + is found to dominate in events with a sin- 
gle light fragment, while Na + dominates in events with a 
multiplicity p of light fragments p > 2. We conclude that 
in the latter events one is observing multifragmentation 
of parent clusters with X > 1. 

We note that the preference for Na + emission for 
X 3> 1 is qualitatively consistent with simple considera- 
tions of the total energy released (Q-value) in a charged 
liquid-drop model which favors the distribution of 
the parent charge over many small fragments. Also, in 
a statistical model that assumes an ergodic distribu- 
tion of fission fragments for a system confined to a small 
volume, evidence was found for a first-order "fragmenta- 
tion phase transition" from a regime with a large residual 
fragment at low excitation energies to a regime with only 
small fragments at high excitation energies. We do not 
observe clear evidence for such an effect in the present 
data, although our statistics and range of study are too 
limited. However, we note that, in addition to providing 
detailed spatial and temporal information, the present 
dynamical approach avoids the ergodic assumption in the 
statistical model. 

Our test system Na24^ + is perhaps too small to be 
able to observe jets unambiguously. We observe only 
small fragments, which tend to be emitted isotropically, 
although the cluster has a tendency to elongate during 
the process (see Figs. || and ||). We stress that these con- 
clusions pertain to the regime of barrierless fission X > 1 . 
We have also searched for symmetric fission when X < 1 , 
where experimental evidence exists for the emission of 



large fission fragments jg]. We found symmetric fission 
for Nai 2+ , and as a rare event (< 10%) for Nais 2+ , 
which has a favorable closed-shell final state Nag + . How- 
ever, in a low-statistics study of Na4o 4+ (A = 1), we ob- 
served only sequential Na 3 + emission (and for Na4o 8+ we 
observed predominantly Na + emission). 
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